Marginal effectsとその周辺

多変量回帰・・・好きですか?

Nozomi Niimi

東京医療センター総合内科

2025-12-30

回帰分析について

皆さん多変量回帰は好きですか?

  • 古くから研究されつくされており、信頼感がある
  • 多くの統計ソフトに入っており、行うのが簡単
  • 解釈性が高く、分かりやすい
  • 本当?

多変量回帰は簡単?

Many regression species

  • 多変量回帰は多くの種類がある
    • 0/1でLogistic回帰
    • 整数値だとPoisson回帰
    • 順序ロジット
    • Censored regression(Cox回帰もこのうち)
  • 選択肢が多く、その分どうすればよいのか分からない!

例えば・・・・・・

  • ICU入室患者に対してRHCが予後を改善するかをみた観察研究(Connors et al. 1996)
    • Propensity scoreを広めた研究としても有名
  • 例えば、RHCが半年以内の死亡と関連するかをロジスティック回帰で解析をする
Parameter Base model
swang1 (RHC) 1.16 [1.01, 1.32]
cat chf (Others) 1.71 [1.29, 2.25]
age 1.03 [1.03, 1.04]
Observations 5733
  • 結果として、RHCは半年以内の死亡と関連する!
  • しかし、ここで上司からのつっこみ

上司からのつっこみ①

結果①

Parameter Base model Interact model
swang1 (RHC) 1.16 [1.01, 1.32] 1.52 [1.02, 2.26]
cat chf (Others) 1.71 [1.29, 2.25] 1.94 [1.40, 2.69]
age 1.03 [1.03, 1.04] 1.03 [1.03, 1.04]
swang1 (RHC) × cat chf (Others) 0.74 [0.49, 1.12]
Observations 5733 5733

上司からのつっこみ②

結果②

Parameter Base model Interact model Spline model
swang1 (RHC) 1.16 [1.01, 1.32] 1.52 [1.02, 2.26] 1.49 [1.00, 2.23]
cat chf (Others) 1.71 [1.29, 2.25] 1.94 [1.40, 2.69] 1.94 [1.40, 2.70]
age 1.03 [1.03, 1.04] 1.03 [1.03, 1.04]
swang1 (RHC) × cat chf (Others) 0.74 [0.49, 1.12] 0.75 [0.49, 1.14]
rcs(age ( degree) 0.99 [0.96, 1.01]
rcs(age ( degree) 1.04 [0.89, 1.21]
rcs(age ( degree) 1.04 [1.03, 1.05]
Observations 5733 5733 5733

最終的に・・・・・・

例えば、2つの連続値をSplineで表したの場合

  • 解釈・・・・・・

我々はどこにいる?

こんな時にJAMA!

我々(臨床医)のしたい事

  • 根本的な問いかけ
    • 「誰への治療効果ですか?」
    • 「どうやった時の治療効果ですか?」

線形回帰の係数だと駄目な理由は?

  • 線形回帰の係数の解釈 = 条件付き期待値
    • 他の値を固定した時に、その因子を1単位変化した時の平均的なアウトカムの変化量
  • その結果、係数(OR, HRなど)という1つの値にまるめてしまう
  • 複雑な関係性をこの係数のみを報告する事で逆にわかりにくくなる事もある
  • Ex. HFmrEFとHFrEFでARNIの効果が変わる

じゃあ、どうやって報告する?

  • 「誰に?」「何が?」「どのくらい?」変化したら、アウトカムが変化するか?
  • 難しい言葉でいうと“Estimand”を報告する必要がある
  • 例えば
    1. 群全体: ATE
    2. 治療を受けた群全体: ATT
    3. 治療を受けなかった群全体: ATU

ここまでの纏め

  • 多変量解析は解釈がわかりにくい!
    • 特に、InteractionやSplineが入るとよりわかりにくい
    • 通常の解析だと、結果は集団全体の平均で丸め込まれてしまう
      • Estimandをどうやって出せばいい?

Standardizationあるいはprametric g-formulaの基礎

Marginal effectsという選択肢

  • G-computationを用いて限界効果を出す手法
  • 本来は、結果のStandardizationの手法
  • Estimandを決定する方法でもある(Hernan 2024)

G-computationの考え方

実践①

  • 生データ
# A tibble: 5,733 × 6
   rowid swang1 estimate conf.low conf.high death_01
   <int> <chr>     <dbl>    <dbl>     <dbl>    <dbl>
 1     1 No RHC    0.606    0.535     0.672        0
 2     2 RHC       0.839    0.787     0.880        1
 3     3 RHC       0.743    0.675     0.801        0
 4     4 No RHC    0.746    0.684     0.800        1
 5     5 RHC       0.872    0.830     0.904        1
 6     6 No RHC    0.780    0.721     0.829        0
 7     7 No RHC    0.585    0.519     0.648        0
 8     8 No RHC    0.287    0.231     0.350        1
 9     9 No RHC    0.315    0.247     0.394        0
10    10 RHC       0.593    0.535     0.648        0
# ℹ 5,723 more rows
  • 元データ+swang1列のみ変更して複製したデータセット
# A tibble: 11,466 × 6
   rowid swang1 estimate conf.low conf.high death_01
   <int> <chr>     <dbl>    <dbl>     <dbl>    <dbl>
 1     1 No RHC    0.606    0.535     0.672        0
 2     2 No RHC    0.823    0.768     0.868        1
 3     3 No RHC    0.721    0.651     0.782        0
 4     4 No RHC    0.746    0.684     0.800        1
 5     5 No RHC    0.859    0.815     0.893        1
 6     6 No RHC    0.780    0.721     0.829        0
 7     7 No RHC    0.585    0.519     0.648        0
 8     8 No RHC    0.287    0.231     0.350        1
 9     9 No RHC    0.315    0.247     0.394        0
10    10 No RHC    0.567    0.508     0.623        0
# ℹ 11,456 more rows
  • 行数が2倍(5733→11466)になっている事に注意!

実践②


 swang1 Estimate Std. Error    z Pr(>|z|)   S 2.5 % 97.5 %
 No RHC    0.638    0.00792 80.5   <0.001 Inf 0.623  0.654
 RHC       0.666    0.01034 64.4   <0.001 Inf 0.645  0.686

Type: response
  • swang1毎でのestimateを纏める
  • 差をとったらRisk differenceを出せる
  • この結果を変形する事で色々な数値を出せる。

実践③

Risk ratio


 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.04   0.0433 4.5     1   1.09

Term: swang1
Type: response
Comparison: ln(mean(RHC) / mean(No RHC))

Odds ratio


 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.13   0.0454 4.5     1   1.27

Term: swang1
Type: response
Comparison: ln(odds(RHC) / odds(No RHC))
  • 疫学的にはリスク差(あるいはリスク比)が最も知りたい数値
  • Odds ratioは発生率が低いイベントでないと、Risk ratioと近似出来ない(Cummings 2009)

G-computtionの応用

G-computationの応用~ATE/ATT/ATU

  • 元データのうち、元々Interventionが0/1の群だけで同様の事をするとATT/ATUも推定可能
  • Interventionだけでなくても、興味がある変数を動かす事で周辺効果(marginal effect)を出すことが可能

実践④-1

  • 例えば、患者背景がCHFでRHCの効果がどうなるかをみたい時
  • まずは交互作用項を含んだ式を作成する
glm(formula = death_01 ~ swang1 + cat_chf + swang1:cat_chf + 
    rcs(age, 4) + crea1 + sex + race + edu + income + wtkilo1 + 
    temp1 + meanbp1 + resp1 + hrt1 + pafi1 + paco21 + ph1 + wblc1 + 
    hema1 + sod1 + pot1 + bili1 + alb1 + cardiohx + chfhx + immunhx + 
    transhx + amihx, family = binomial, data = rhc_prep)
Parameter Coefficient SE 95% CI z p
(Intercept) 754.71 2013.26 (4.14, 1.45e+05) 2.48 0.013
swang1 (RHC) 1.49 0.30 (1.00, 2.23) 1.97 0.049
cat chf (Others) 1.94 0.33 (1.40, 2.70) 3.97 < .001
rcs(age ( degree) 1.04 6.43e-03 (1.03, 1.05) 6.24 < .001
rcs(age ( degree) 0.99 0.01 (0.96, 1.01) -0.98 0.325
rcs(age ( degree) 1.04 0.08 (0.89, 1.21) 0.51 0.611
crea1 1.02 0.02 (0.99, 1.05) 1.16 0.244
sex (Male) 1.23 0.08 (1.09, 1.39) 3.42 < .001
race (other) 1.05 0.15 (0.80, 1.38) 0.33 0.741
race (white) 0.95 0.08 (0.80, 1.12) -0.61 0.542
edu 1.02 0.01 (1.00, 1.04) 1.66 0.098
income ($11-$25k) 1.15 0.14 (0.90, 1.47) 1.13 0.257
income ($25-$50k) 1.00 0.13 (0.78, 1.28) 0.02 0.983
income (Under $11k) 1.40 0.17 (1.10, 1.77) 2.77 0.006
wtkilo1 1.00 1.05e-03 (0.99, 1.00) -3.80 < .001
temp1 0.92 0.02 (0.89, 0.95) -4.49 < .001
meanbp1 1.00 8.02e-04 (0.99, 1.00) -4.73 < .001
resp1 1.00 2.20e-03 (1.00, 1.01) 1.84 0.066
hrt1 1.00 7.87e-04 (1.00, 1.00) 3.23 0.001
pafi1 1.00 2.82e-04 (1.00, 1.00) 1.85 0.065
paco21 0.99 2.74e-03 (0.99, 1.00) -2.95 0.003
ph1 0.46 0.16 (0.23, 0.89) -2.29 0.022
wblc1 1.00 2.59e-03 (1.00, 1.01) 0.58 0.564
hema1 0.98 3.94e-03 (0.97, 0.99) -4.82 < .001
sod1 1.00 3.99e-03 (0.99, 1.01) 0.50 0.617
pot1 1.03 0.03 (0.96, 1.09) 0.79 0.431
bili1 1.06 9.33e-03 (1.04, 1.08) 6.43 < .001
alb1 0.96 0.04 (0.88, 1.05) -0.82 0.411
cardiohx 1.08 0.10 (0.90, 1.29) 0.83 0.406
chfhx 1.52 0.16 (1.25, 1.86) 4.13 < .001
immunhx 1.25 0.09 (1.09, 1.44) 3.26 0.001
transhx 0.70 0.06 (0.59, 0.84) -3.94 < .001
amihx 0.78 0.12 (0.57, 1.07) -1.55 0.122
swang1 (RHC) × cat chf (Others) 0.75 0.16 (0.49, 1.13) -1.37 0.172

実践④-2

  • 次に、データセットを介入に併せて(RHCの有無)拡張した上でSubgroupを作成する

実践④-3

  • このsubgroupにモデルをfitさせて個々人の死亡率を計算する
# A tibble: 912 × 9
   rowid rowidcf swang1 cat_chf   age crea1 estimate conf.low conf.high
   <int>   <int> <chr>  <chr>   <dbl> <dbl>    <dbl>    <dbl>     <dbl>
 1     1       1 No RHC CHF      42.1 1.30     0.360    0.266     0.467
 2     2       2 No RHC CHF      62.1 1.30     0.610    0.534     0.681
 3     3       3 No RHC CHF      73.5 1.90     0.687    0.614     0.752
 4     4       4 No RHC CHF      41.3 1.5      0.335    0.260     0.420
 5     5       5 No RHC CHF      58.0 1.60     0.457    0.361     0.557
 6     6       6 No RHC CHF      68.4 0.800    0.475    0.387     0.564
 7     7       7 No RHC CHF      63.5 1.10     0.561    0.472     0.645
 8     8       8 No RHC CHF      66.9 2.40     0.566    0.476     0.652
 9     9       9 No RHC CHF      57.5 1.10     0.597    0.514     0.675
10    10      10 No RHC CHF      76.4 2        0.673    0.568     0.764
# ℹ 902 more rows

 - 例えば、平均するとこんな感じ 

 swang1 Estimate Std. Error    z Pr(>|z|)     S 2.5 % 97.5 %
 No RHC    0.562     0.0306 18.3   <0.001 247.4 0.502  0.622
 RHC       0.651     0.0317 20.5   <0.001 307.8 0.588  0.713

Type: response

実践④-4

  • このEstimateを使っていろんな指標を算出する

    • Risk differenceの場合

 Estimate Std. Error    z Pr(>|z|)   S   2.5 % 97.5 %
   0.0884     0.0444 1.99   0.0466 4.4 0.00132  0.175

Term: swang1
Type: response
Comparison: RHC - No RHC
  • Odds ratioの場合

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.16   0.0474 4.4     1   1.34

Term: swang1
Type: response
Comparison: ln(mean(RHC) / mean(No RHC))

G-computationの応用~Propensity score match

  1. Propensity score matchを計算
  2. Matched cohortを作成
  3. その上で、共変量を入れたアウトカムモデルで回帰を行う
  4. RCTの多変量と同じでバランスをちゃんと取れる
  5. どの群を選ぶかでATT/ATE/ATUも簡単に計算可能!(Greifer 2025)

実践⑤-1

  • まずはPropensity score match行い、Matched dataを出す
A `matchit` object
 - method: 1:1 nearest neighbor matching without replacement
 - distance: Propensity score
             - estimated with logistic regression
 - number of obs.: 5733 (original), 4366 (matched)
 - target estimand: ATT
 - covariates: cat_chf, age, sex, race, edu, income, wtkilo1, temp1, meanbp1, resp1, hrt1, pafi1, paco21, ph1, wblc1, hema1, sod1, pot1, crea1, bili1, alb1, cardiohx, chfhx, immunhx, transhx, amihx
matched data 
# A tibble: 4,366 × 6
   death_yn swang1 cat_chf   age distance weights
      <dbl> <chr>  <chr>   <dbl>    <dbl>   <dbl>
 1        0 No RHC Others   70.3    0.502       1
 2        1 RHC    Others   78.2    0.563       1
 3        0 RHC    Others   46.1    0.402       1
 4        1 No RHC Others   75.3    0.344       1
 5        1 RHC    Others   67.9    0.302       1
 6        0 No RHC Others   55.0    0.379       1
 7        1 No RHC Others   43.6    0.281       1
 8        0 No RHC Others   18.0    0.283       1
 9        0 RHC    Others   48.4    0.490       1
10        0 No RHC Others   34.4    0.393       1
# ℹ 4,356 more rows

実践⑤-2

  • Propensity scoreに用いた共変量を含めてアウトカムの回帰モデルを作成する
  • この重みも用いることに注意!
アウトカムモデル 
glm(formula = death_01 ~ swang1 * cat_chf + rcs(age, 4) + crea1 + 
    sex + race + edu + income + wtkilo1 + temp1 + meanbp1 + resp1 + 
    hrt1 + pafi1 + paco21 + ph1 + wblc1 + hema1 + sod1 + pot1 + 
    bili1 + alb1 + cardiohx + chfhx + immunhx + transhx + amihx, 
    family = binomial, data = dat_m, weights = weights)

実践⑤~結果

  • G-computationを各々の群で行う
# A tibble: 4,366 × 8
   rowid contrast     estimate conf.low conf.high death_01 swang1 weights
   <int> <chr>           <dbl>    <dbl>     <dbl>    <dbl> <chr>    <dbl>
 1     1 RHC - No RHC  0.0196  -0.0137     0.0529        0 No RHC       1
 2     2 RHC - No RHC  0.0120  -0.00850    0.0324        1 RHC          1
 3     3 RHC - No RHC  0.0163  -0.0115     0.0441        0 RHC          1
 4     4 RHC - No RHC  0.0156  -0.0110     0.0422        1 No RHC       1
 5     5 RHC - No RHC  0.00965 -0.00700    0.0263        1 RHC          1
 6     6 RHC - No RHC  0.0193  -0.0134     0.0519        0 No RHC       1
 7     7 RHC - No RHC  0.0169  -0.0120     0.0458        1 No RHC       1
 8     8 RHC - No RHC  0.0177  -0.0126     0.0480        0 No RHC       1
 9     9 RHC - No RHC  0.0201  -0.0140     0.0543        0 RHC          1
10    10 RHC - No RHC  0.0197  -0.0139     0.0533        0 No RHC       1
# ℹ 4,356 more rows


ここから計算する

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.03    0.122 3.0 0.992   1.07

Term: swang1
Type: response
Comparison: ln(mean(RHC) / mean(No RHC))
# A tibble: 2,183 × 8
   rowid contrast     estimate conf.low conf.high death_01 swang1 weights
   <int> <chr>           <dbl>    <dbl>     <dbl>    <dbl> <chr>    <dbl>
 1     1 RHC - No RHC  0.0120  -0.00850    0.0324        1 RHC          1
 2     2 RHC - No RHC  0.0163  -0.0115     0.0441        0 RHC          1
 3     3 RHC - No RHC  0.00965 -0.00700    0.0263        1 RHC          1
 4     4 RHC - No RHC  0.0201  -0.0140     0.0543        0 RHC          1
 5     5 RHC - No RHC  0.0207  -0.0145     0.0560        0 RHC          1
 6     6 RHC - No RHC  0.0199  -0.0140     0.0538        1 RHC          1
 7     7 RHC - No RHC  0.00868 -0.00604    0.0234        1 RHC          1
 8     8 RHC - No RHC  0.0165  -0.0115     0.0446        1 RHC          1
 9     9 RHC - No RHC  0.0146  -0.0106     0.0399        0 RHC          1
10    10 RHC - No RHC  0.0206  -0.0143     0.0556        1 RHC          1
# ℹ 2,173 more rows


ここから計算する

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.03    0.121 3.0 0.992   1.07

Term: swang1
Type: response
Comparison: ln(mean(RHC) / mean(No RHC))
# A tibble: 2,183 × 8
   rowid contrast     estimate conf.low conf.high death_01 swang1 weights
   <int> <chr>           <dbl>    <dbl>     <dbl>    <dbl> <chr>    <dbl>
 1     1 RHC - No RHC   0.0196 -0.0137     0.0529        0 No RHC       1
 2     2 RHC - No RHC   0.0156 -0.0110     0.0422        1 No RHC       1
 3     3 RHC - No RHC   0.0193 -0.0134     0.0519        0 No RHC       1
 4     4 RHC - No RHC   0.0169 -0.0120     0.0458        1 No RHC       1
 5     5 RHC - No RHC   0.0177 -0.0126     0.0480        0 No RHC       1
 6     6 RHC - No RHC   0.0197 -0.0139     0.0533        0 No RHC       1
 7     7 RHC - No RHC   0.0210 -0.0146     0.0565        1 No RHC       1
 8     8 RHC - No RHC   0.0203 -0.0143     0.0550        0 No RHC       1
 9     9 RHC - No RHC   0.0124 -0.00937    0.0342        0 No RHC       1
10    10 RHC - No RHC   0.0642 -0.0284     0.157         1 No RHC       1
# ℹ 2,173 more rows
ここから計算する

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.03    0.123 3.0 0.992   1.07

Term: swang1
Type: response
Comparison: ln(mean(RHC) / mean(No RHC))

G-computationの応用~Doubly robust

  1. Propensity score weightingを計算
  2. Outcomeを目標とする多変量回帰を作成
  3. 上記を組み合わせてDoubly robustを計算可能
  4. どの群を選ぶかでATT/ATE/ATUも簡単に計算可能!

実践⑤~結果

  1. treatmentを目的としたLogistic regressionを作成して、予測式を作成
  2. 予測式とATE/ATT/ATCでそれぞれ重み付けの計算方法を変える(ATOなどもあり)
  3. この重みを考慮したアウトカムモデルを作成する
  4. このアウトカムモデルで、各データでのG-computationを行う
  5. ATEは全体、ATTはRHC使用者、ATCはRHC不使用者
 
 1. 重み付けの式 
WeightIt::weightit(formula = swang_yn ~ cat_chf + age + sex + 
    race + edu + income + wtkilo1 + temp1 + meanbp1 + resp1 + 
    hrt1 + pafi1 + paco21 + ph1 + wblc1 + hema1 + sod1 + pot1 + 
    crea1 + bili1 + alb1 + cardiohx + chfhx + immunhx + transhx + 
    amihx, data = rhc_prep, method = "glm", estimand = "ATE")

 2. 重みを元のデータセットに追加 
# A tibble: 5,733 × 8
   death_01 swang_yn   age sex    race  cat_chf crea1 weights
      <dbl>    <dbl> <dbl> <chr>  <chr> <chr>   <dbl>   <dbl>
 1        0        0  70.3 Male   white Others  1.20     2.01
 2        1        1  78.2 Female white Others  0.600    1.78
 3        0        1  46.1 Female white Others  2.60     2.49
 4        1        0  75.3 Female white Others  1.70     1.53
 5        1        1  67.9 Male   white Others  3.60     3.31
 6        0        0  86.1 Female white Others  1.40     1.12
 7        0        0  55.0 Male   white Others  1        1.61
 8        1        0  43.6 Male   white Others  0.700    1.39
 9        0        0  18.0 Female white Others  1.70     1.39
10        0        1  48.4 Female white Others  0.5      2.04
# ℹ 5,723 more rows

 3. 重みを考慮したアウトカムモデル 
WeightIt::glm_weightit(formula = death_01 ~ swang1 * cat_chf + 
    rcs(age, 4) + crea1 + sex + race + edu + income + wtkilo1 + 
    temp1 + meanbp1 + resp1 + hrt1 + pafi1 + paco21 + ph1 + wblc1 + 
    hema1 + sod1 + pot1 + bili1 + alb1 + cardiohx + chfhx + immunhx + 
    transhx + amihx, data = rhc_prep, family = binomial, weightit = wout_ate)

 4. G-computation 
# A tibble: 11,466 × 6
   rowid   age swang1 estimate conf.low conf.high
   <int> <dbl> <chr>     <dbl>    <dbl>     <dbl>
 1     1  70.3 No RHC    0.562    0.474     0.647
 2     2  78.2 No RHC    0.801    0.729     0.858
 3     3  46.1 No RHC    0.684    0.590     0.765
 4     4  75.3 No RHC    0.756    0.678     0.820
 5     5  67.9 No RHC    0.864    0.813     0.903
 6     6  86.1 No RHC    0.824    0.757     0.876
 7     7  55.0 No RHC    0.566    0.482     0.646
 8     8  43.6 No RHC    0.258    0.200     0.325
 9     9  18.0 No RHC    0.343    0.251     0.449
10    10  48.4 No RHC    0.547    0.476     0.616
# ℹ 11,456 more rows

 5. Risk ratio and 95% confidence interval 

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.05   0.0291 5.1     1    1.1

Term: swang1
Type: probs
Comparison: ln(mean(RHC) / mean(No RHC))
 
 1. 重み付けの式 
WeightIt::weightit(formula = swang_yn ~ cat_chf + age + sex + 
    race + edu + income + wtkilo1 + temp1 + meanbp1 + resp1 + 
    hrt1 + pafi1 + paco21 + ph1 + wblc1 + hema1 + sod1 + pot1 + 
    crea1 + bili1 + alb1 + cardiohx + chfhx + immunhx + transhx + 
    amihx, data = rhc_prep, method = "glm", estimand = "ATT")

 2. 重みを元のデータセットに追加 
# A tibble: 5,733 × 8
   death_01 swang_yn   age sex    race  cat_chf crea1 weights
      <dbl>    <dbl> <dbl> <chr>  <chr> <chr>   <dbl>   <dbl>
 1        0        0  70.3 Male   white Others  1.20    1.01 
 2        1        1  78.2 Female white Others  0.600   1    
 3        0        1  46.1 Female white Others  2.60    1    
 4        1        0  75.3 Female white Others  1.70    0.526
 5        1        1  67.9 Male   white Others  3.60    1    
 6        0        0  86.1 Female white Others  1.40    0.117
 7        0        0  55.0 Male   white Others  1       0.611
 8        1        0  43.6 Male   white Others  0.700   0.392
 9        0        0  18.0 Female white Others  1.70    0.394
10        0        1  48.4 Female white Others  0.5     1    
# ℹ 5,723 more rows

 3. 重みを考慮したアウトカムモデル 
WeightIt::glm_weightit(formula = death_01 ~ swang1 * cat_chf + 
    rcs(age, 4) + crea1 + sex + race + edu + income + wtkilo1 + 
    temp1 + meanbp1 + resp1 + hrt1 + pafi1 + paco21 + ph1 + wblc1 + 
    hema1 + sod1 + pot1 + bili1 + alb1 + cardiohx + chfhx + immunhx + 
    transhx + amihx, data = rhc_prep, family = binomial, weightit = wout_att)

 4. G-computation 
# A tibble: 4,366 × 6
   rowid   age swang1 estimate conf.low conf.high
   <int> <dbl> <chr>     <dbl>    <dbl>     <dbl>
 1     1  78.2 No RHC    0.834    0.770     0.883
 2     2  46.1 No RHC    0.733    0.643     0.808
 3     3  67.9 No RHC    0.864    0.810     0.905
 4     4  48.4 No RHC    0.582    0.513     0.649
 5     5  68.3 No RHC    0.554    0.434     0.669
 6     6  74.7 No RHC    0.604    0.506     0.694
 7     7  88.4 No RHC    0.881    0.825     0.920
 8     8  69.0 No RHC    0.719    0.656     0.775
 9     9  50.6 No RHC    0.800    0.716     0.864
10    10  62.7 No RHC    0.577    0.493     0.657
# ℹ 4,356 more rows

 5. Risk ratio and 95% confidence interval 

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.04    0.115 3.1 0.991   1.08

Term: swang1
Type: probs
Comparison: ln(mean(RHC) / mean(No RHC))
 
 1. 重み付けの式 
WeightIt::weightit(formula = swang_yn ~ cat_chf + age + sex + 
    race + edu + income + wtkilo1 + temp1 + meanbp1 + resp1 + 
    hrt1 + pafi1 + paco21 + ph1 + wblc1 + hema1 + sod1 + pot1 + 
    crea1 + bili1 + alb1 + cardiohx + chfhx + immunhx + transhx + 
    amihx, data = rhc_prep, method = "glm", estimand = "ATC")

 2. 重みを元のデータセットに追加 
# A tibble: 5,733 × 8
   death_01 swang_yn   age sex    race  cat_chf crea1 weights
      <dbl>    <dbl> <dbl> <chr>  <chr> <chr>   <dbl>   <dbl>
 1        0        0  70.3 Male   white Others  1.20    1    
 2        1        1  78.2 Female white Others  0.600   0.777
 3        0        1  46.1 Female white Others  2.60    1.49 
 4        1        0  75.3 Female white Others  1.70    1    
 5        1        1  67.9 Male   white Others  3.60    2.31 
 6        0        0  86.1 Female white Others  1.40    1    
 7        0        0  55.0 Male   white Others  1       1    
 8        1        0  43.6 Male   white Others  0.700   1    
 9        0        0  18.0 Female white Others  1.70    1    
10        0        1  48.4 Female white Others  0.5     1.04 
# ℹ 5,723 more rows

 3. 重みを考慮したアウトカムモデル 
WeightIt::glm_weightit(formula = death_01 ~ swang1 * cat_chf + 
    rcs(age, 4) + crea1 + sex + race + edu + income + wtkilo1 + 
    temp1 + meanbp1 + resp1 + hrt1 + pafi1 + paco21 + ph1 + wblc1 + 
    hema1 + sod1 + pot1 + bili1 + alb1 + cardiohx + chfhx + immunhx + 
    transhx + amihx, data = rhc_prep, family = binomial, weightit = wout_atc)

 4. G-computation 
# A tibble: 7,100 × 6
   rowid   age swang1 estimate conf.low conf.high
   <int> <dbl> <chr>     <dbl>    <dbl>     <dbl>
 1     1  70.3 No RHC    0.553    0.454     0.648
 2     2  75.3 No RHC    0.771    0.685     0.839
 3     3  86.1 No RHC    0.828    0.753     0.884
 4     4  55.0 No RHC    0.530    0.435     0.623
 5     5  43.6 No RHC    0.240    0.179     0.313
 6     6  18.0 No RHC    0.401    0.288     0.527
 7     7  34.4 No RHC    0.347    0.258     0.448
 8     8  42.2 No RHC    0.436    0.330     0.547
 9     9  82.0 No RHC    0.630    0.513     0.733
10    10  78.3 No RHC    0.747    0.673     0.808
# ℹ 7,090 more rows

 5. Risk ratio and 95% confidence interval 

 Estimate Pr(>|z|)   S 2.5 % 97.5 %
     1.06   0.0271 5.2  1.01   1.12

Term: swang1
Type: probs
Comparison: ln(mean(RHC) / mean(No RHC))

G-computationの応用~機械学習を用いた手法

  • 厳密にはStandardizationとは違うが・・・・・・
  • Datasetを倍にして機械学習で推定する方法をS-learnerという
  • 因果推定を行う場合は、重み付けも行うD-learnerやX-learnerもある
  • 信頼区間やP値が出せないのが難点

実践⑤ S-learnerの場合

実践⑤~モデル作成

  • Machine learning modelとしてXGBoostとする
  • ハイパーパラメーターはデフォルト
  • 何も考えずに、全体で学習させる
モデルの中身
Boosted Tree Model Specification (classification)

Computational engine: xgboost 
##### xgb.Booster
call:
  xgboost::xgb.train(params = list(eta = 0.3, max_depth = 6, gamma = 0, 
    colsample_bytree = 1, colsample_bynode = 1, min_child_weight = 1, 
    subsample = 1, nthread = 1, objective = "binary:logistic"), 
    data = x$data, nrounds = 15, evals = x$watchlist, verbose = 0)
# of features: 38 
# of rounds:  15 
callbacks:
   evaluation_log 
evaluation_log:
  iter training_logloss
 <int>            <num>
     1        0.6098085
     2        0.5866390
   ---              ---
    14        0.4528868
    15        0.4455501

実践⑤ 仮想データ作成

# A tibble: 11,466 × 33
   rowid rowidcf death_yn death_days swang_yn cat_chf cat1       age crea1 sex  
   <int>   <int>    <dbl>      <dbl>    <dbl> <chr>   <chr>    <dbl> <dbl> <chr>
 1     1       1        0        180        0 Others  COPD      70.3 1.20  Male 
 2     2       2        1         45        1 Others  MOSF w/…  78.2 0.600 Fema…
 3     3       3        0        180        1 Others  MOSF w/…  46.1 2.60  Fema…
 4     4       4        1         37        0 Others  ARF       75.3 1.70  Fema…
 5     5       5        1          2        1 Others  MOSF w/…  67.9 3.60  Male 
 6     6       6        0        180        0 Others  COPD      86.1 1.40  Fema…
 7     7       7        0        180        0 Others  MOSF w/…  55.0 1     Male 
 8     8       8        1         38        0 Others  ARF       43.6 0.700 Male 
 9     9       9        0        180        0 Others  MOSF w/…  18.0 1.70  Fema…
10    10      10        0        180        1 Others  ARF       48.4 0.5   Fema…
# ℹ 11,456 more rows
# ℹ 23 more variables: race <chr>, edu <dbl>, income <chr>, wtkilo1 <dbl>,
#   temp1 <dbl>, meanbp1 <dbl>, resp1 <dbl>, hrt1 <int>, pafi1 <dbl>,
#   paco21 <dbl>, ph1 <dbl>, wblc1 <dbl>, hema1 <dbl>, sod1 <int>, pot1 <dbl>,
#   bili1 <dbl>, alb1 <dbl>, cardiohx <int>, chfhx <int>, immunhx <int>,
#   transhx <int>, amihx <int>, swang1 <chr>

実践⑤ 結果を纏める

# A tibble: 11,466 × 4
   swang1 death_01 probability_1 predicted_class
   <chr>  <fct>            <dbl> <fct>          
 1 No RHC 0                0.601 1              
 2 No RHC 1                0.795 1              
 3 No RHC 0                0.545 1              
 4 No RHC 1                0.745 1              
 5 No RHC 1                0.859 1              
 6 No RHC 0                0.642 1              
 7 No RHC 0                0.599 1              
 8 No RHC 1                0.265 0              
 9 No RHC 0                0.546 1              
10 No RHC 0                0.426 0              
# ℹ 11,456 more rows
# A tibble: 2 × 2
  swang1 mean_death_prob
  <chr>            <dbl>
1 No RHC           0.644
2 RHC              0.654

G-computationの利点

  • 「この集団の介入を変えたら、どの程度良くなるか?」をダイレクトに伝えられる(King, Tomz, and Wittenberg 2000)
    • InteractionやSplineなど複雑な式でもシンプルに結果を伝えられる
  • 予測と因果を両方行う事が可能!

予測~Average marginal prediction

  • 通常のアウトカム式のみで一発勝負
  • 一応Bias補正方法も存在する(Tibshirani et al. 2024)

因果~Average comparison

  • AIPWやTLMEを使った方が安心かもしれない
    • ただし、効率が悪い可能性もあり

G-computationでのAdvancedな統計結果の伝え方

  • 伝えるEstimandを正確に伝える
  • 「誰に」、「もし**したら」、「平均的にどのような効果?」が「どれくらいの確実性」あるかを言う
    • 出来るだけ、臨床的に意味がある差かどうかを考える

G-computationの限界

  • 基本的には通常の回帰と一緒
    • Unobserved confoundingやMisspecificationに弱い
  • 予測の為ではなく因果よりと考えたほうが良い
    • モデル作成でStep wiseや機械学習のようなやり方はしないほうが良い

自分の研究においては・・・・・・

  • AF患者において動悸およびそれ以外の症状(息切れ、めまい/失神、倦怠感)のHR-QOLへの影響をみた研究
  • 特に、非動悸の患者においてもCAは有効なのか?を疑問として提案
  • 動悸、息切れ、めまい/失神、倦怠感の有無でCAの有効性をみたかった

How to use marginal effects

  • これを、動悸/めまい/息切れ/倦怠感で各々行った
    • さらに、感度分析で機械学習モデル(grf)を使ってみたりしてみた

Result

forest plot

responder fig

ML fig
Figure 1: Main figure
  • 誰に(症状がある患者)、どの介入(Rate/AAD/CA)を行うとQOLがどうなるか?を臨床的に使いやすく見せることを目的
  • そのため、Responder analysisも行った
  • 機械学習の結果もRobustnessを示すために行った

最後に

Box先生の名言

すべてのモデルは誤っている。しかし、そのうちのいくつかは役に立つ。

RでのMarginal effectsの使い方

  • emmeans, marginaleffects, easystatsmodelbasedなど
  • AIMの論文にも記載があり使いやすいのはmarginaleffects

Thank you for your listening!!

References

Connors, A F, Jr, T Speroff, N V Dawson, C Thomas, F E Harrell Jr, D Wagner, N Desbiens, et al. 1996. “The Effectiveness of Right Heart Catheterization in the Initial Care of Critically Ill Patients. SUPPORT Investigators.” JAMA: The Journal of the American Medical Association 276 (11): 889–97. https://doi.org/10.1001/jama.276.11.889.
Cummings, Peter. 2009. “The Relative Merits of Risk Ratios and Odds Ratios.” Archives of Pediatrics & Adolescent Medicine 163 (5): 438–45. https://doi.org/10.1001/archpediatrics.2009.31.
Dahabreh, Issa J, and Kirsten Bibbins-Domingo. 2024. “Causal Inference about the Effects of Interventions from Observational Studies in Medical Journals.” JAMA: The Journal of the American Medical Association 331 (21): 1845–53. https://doi.org/10.1001/jama.2024.7741.
Greifer, Noah. 2025. “Estimating Effects After Matching.” https://cran.r-project.org/web/packages/MatchIt/vignettes/estimating-effects.html.
Hernan, Miguel A. 2024. Causal Inference: What If. Edited by James M. Robins. Boca Raton: Taylor & Francis.
King, Gary, Michael Tomz, and Jason Wittenberg. 2000. “Making the Most of Statistical Analyses: Improving Interpretation and Presentation.” American Journal of Political Science 44 (2): 347–61. https://doi.org/10.2307/2669316.
Tibshirani, Julie, Susan Athey, Erik Sverdrup, and Stefan Wager. 2024. Grf: Generalized Random Forests. https://github.com/grf-labs/grf.